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THEORY  OF  VORTICITY  GENERATION  BY  SHOCK  WAVE 
AND  FLAME  INTERACTIONS 


INTRODUCTION 

The  interaction  of  pressure  waves  with  density  gradients  is  a  fundamental 
source  of  long-lived  vorticity  in  fluids  (Chu  and  Kovasznay,  1957,  and  Picone 
and  Boris,  1983).  This  mechanism  is  particularly  important  in  combustion, 
since  the  release  of  chemical  energy  produces  both  pressure  and  density 
disturbances  in  the  fluid.  These  disturbances  then  interact,  producing 
significant  vorticity.  The  rotational  motion  associated  with  the  vorticity 
field  enhances  mixing  and  introduces  additional  fluctuations  in  the  flow 
variables.  The  scale  lengths  which  characterize  these  secondary  fluctuations 
are  roughly  one-half  those  of  the  local  density  gradients  which  exist  prior  to 
interacting  with  given  pressure  waves  (Picone  and  Boris,  1983).  Because  the 
structures  and  sizes  of  the  fluid  interfaces  change  under  the  influence  of  the 
newly  generated  vorticity,  local  reaction  rates  change,  amplifying  the  effects 
of  the  fluctuations  and  producing  more  local  pressure  waves.  The  process  thus 
continues,  progressively  reducing  the  length  scales  of  inhomogeneities  in  the 
flow  field.  Through  this  mechanism,  a  flow  which  is  originally  laminar  can 
quickly  become  turbulent.  Such  phenomena  occur  in  almost  every  nonidealized 
high  Reynolds  number  flow. 

In  this  paper,  we  focus  on  the  interaction  of  a  weak,  planar  shock  wave 

with  an  azimuthally  symmetric  region  characterized  by  density  and  temperature 

distributions  similar  to  those  found  in  an  expanding  flame.  Our  calculations 

are  two-dimensional  and  correspond  to  a  cylindrical  flame.  This  problem 

contains  the  essential  features  of  vorticity  generation  by  pressure  wave  and 

flame  interactions.  Ignoring  the  flame  dynamics  is  acceptable  for  calculating 

the  large-scale  vorticity  distribution,  since  the  transit  time  of  the  shock 
Manuscript  approved  February  1,  1984. 
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across  the  heated  region  is  quite  short  compared  to  the  time  scale  of  flame 
propagation.  Markstein  has  studied  a  similar  case  experimentally  for  a 
reactive  medium,  using  a  long,  vertical  shock  tube  with  a  30  cm  combustion 
chamber  at  the  bottom  and  a  diaphragm  90  cm  from  the  top.  Figure  1  shows 
Schlieren  photographs  from  one  of  Markstein' s  experiments.  A  weak  shock 
(pressure  ratio  ~  1.3)  passes  through  a  roughly  spherical  flame  approximately 
15  cm  from  the  bottom  of  the  chamber,  which  contains  a  stoichiometric  mixture 
of  n-butane  and  air.  In  the  first  frame,  we  see  the  incident  shock  less  than 
1  cm  from  the  flame  boundary.  The  flame  actually  appears  to  be  more  oblong 
than  spherical.  In  the  second  photograph,  for  which  time  is  defined  to  be 
0.10  ms,  we  see  a  curved  rarefaction  wave  moving  upward  from  the  flame  while 
the  upper  flame  boundary  has  been  compressed  by  the  flow  behind  the  incident 
shock.  By  0.40  ms,  a  vortex  ring  has  formed  due  to  the  interaction  of  the 
shock  wave  and  the  flame.  The  enhanced  flow  at  the  center  of  the  ring  pulls 
unreacted  gas  through  the  flame,  causing  the  gas  to  ignite.  As  the  newly 
ignited  gas  emerges,  a  fine-grained  turbulent  burning  zone  develops  (apparent 
at  0.70  ms ) . 

To  analyze  the  fluid-dynamic  aspects  of  vorticity  generation,  we  will 
study  the  closely  related  situation  in  which  a  planar  shock  interacts  with  a 
hot,  low  density  region  in  a  nonreactive  fluid.  Past  treatments  of  this 
problem  have  relied  on  the  linearization  of  all  or  part  of  the  relevant 
equations  or  have  used  a  perturbation  expansion  (Rudinger,  1958,  and 
Markstein,  1964).  Those  approaches  are  valid  for  cases  with  vanishingly  small 
perturbation  amplitudes  and  time  scales  but  cannot  provide  a  quantitative 
picture  of  the  large  scale  flows  which  result  from  the  Initial  asymmetries, 
for  which  the  characteristic  length  scale  is  the  flame  radius  and  the 
interaction  time  scale  Is  finite. 
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In  the  next  section,  we  outline  a  simple  nonlinear  theory  of  vorticity 
generation  by  shock  wave  and  flame  interactions.  We  then  present  the  results 
of  a  two-dimensional  numerical  simulation  of  the  nonreactive  case.  The 
simulation  provides  a  clear  picture  of  the  formation  of  a  vortex  ring.  The 
calculation  also  provides  a  calibration  of  constants  arising  in  the  theory. 
Finally,  we  discuss  the  shock  wave  and  flame  interaction  as  a  mechanism  for 
developing  and  enhancing  turbulence  and  for  determining  the  set  of  turbulent 
scales  observed  experimentally. 
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NONLINEAR  THEORY 


Equation  (1)  gives  the  rigorous,  inviscid  equation  for  the  time  evolution 
of  vorticity, 

d£ 

dt’+£~*Y.=  £*£Y.+  (ZpxHP)/p2>  ( i ) 

where  v  is  the  fluid  velocity,  £  =  V  x  v  is  the  vorticity,  p  is  the  density, 
and  P  is  the  pressure.  All  variables  are  functions  of  the  position  r  and 
time  t.  Equation  (1)  provides  a  direct  mechanism  for  the  generation  of 
vorticity  by  shock  wave  and  flame  interactions,  since  the  right  hand  side 
contains  a  source  term  which  is  proportional  to  the  cross  product  of  the 
density  and  pressure  gradients.  Whenever  the  local  pressure  and  density 
gradients  are  misaligned,  the  source  term  will  be  nonzero  and  the  production 
of  vorticity  will  occur.  In  Fig.  1,  the  experimental  case  of  interest,  the 
pressure  gradient  associated  with  the  planar  shock  is  parallel  to  the  axis  of 
the  shock  tube  while  che  density  gradient  of  the  spherical  flame  is  directed 
approximately  radially  outward  from  the  center.  Large  scale  vorticity  in  the 
form  of  a  vortex  ring  results  from  the  interaction. 

An  estimate  of  vorticity  generation  in  the  nonreactive  case  by  this 
mechanism  requires  the  integration  of  Eq.  (1)  in  both  space  and  time.  We 
may,  however,  simplify  our  task  by  working  in  two  dimensions  (Cartesian), 
since  the  induced  rotational  flows,  which  are  cylindrically  symmetric,  differ 
from  those  of  a  vortex  ring  by  geometrical  factors  of  order  one  (see  sentence 
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following  Eq.  (5)  below).  We  then  exploit  the  fact  that  after  the  shock  has 
passed  through  a  circular  region  of  density  variation  in  the  x-y  plane,  the 
residual  flow  field  consists  of  a  pair  of  vortices  of  strength  -  <z  =  -  < 
(Picone  and  Boris,  1983).  Figure  2  defines  the  notation  for  this  section  and 
the  Appendix.  For  the  purpose  of  integrating  Eq.  (1)  analytically,  we  align 
the  x-axis  with  the  direction  of  propagation  of  the  shock  and  place  the  origin 
at  the  center  of  the  density  depression.  The  quantity  Sq  denotes  the  radius 
of  the  heated  region,  and  the  vortices  will  be  centered  at  (x,y)  =>  (0,  ±  y), 
where  y  ^-Sg.  That  is,  in  the  case  of  a  planar  shock  propagating  downward 
through  a  circular  heated  region  in  a  two-dimensional  vertical  chamber,  there 
will  be  one  vortex  on  either  side  of  the  vertical  axis  of  the  chamber.  The 
pair  of  vorticies  will  be  roughly  centered  on  the  heated  region,  and  the 
associated  rotational  flows  /ill  be  oppositely  directed.  The  strength  of  one 
vortex  at  time  t  is  equal  to  the  integral  of  the  vorticity  over  the  half  plane 
containing  the  vortex.  For  the  case  of  interest,  the  vortex  strength  is 

CO  30 

3  l dy  $ dx  sz<x»y»o» 

”  ~  CO 

which  satisfies  the  equation  (see  the  Appendix) 

d  x  00  00 

=»  /  dy  /  dx  [(?p  x  VP )z/  P2]  • 

0  —00 

In  a  formal  sense,  the  residual  vortex  strength  is  then 

k  =  lira  <  ( t ) .  (4) 

t  >  ®  z 


(2) 


(3) 


& 
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In  practice,  the  generation  of  vorticity  occurs  over  a  finite  time 
interval  (-t/2,  t/2).  For  an  azimuthally  symmetric  density  depresssion  with  a 
characteristic  radius  Sq,  the  time  interval  t  is  approximately  equal  to  2S0/W, 
where  W  is  the  average  velocity  with  which  the  shock  propagates  across  the 
region.  Because  this  time  is  much  shorter  than  the  time  required  for  the 
resulting  rotational  flows  to  affect  the  density  gradient  appreciably,  we  may 
ignore  the  feedback  effects  of  the  vorticity  on  the  generation  process. 

Given  the  value  of  k,  we  may  compute  the  fluid  velocity  along  the  line 
bisecting  the  vortex  pair  by  summing  the  flows  induced  by  each  vortex 
separately.  Because  the  azimuthal  velocity  of  an  isolated  vortex  filament  has 
the  magnitude  |vg  |  =  j <  |/ 2 irr ,  where  r  is  the  radial  distance  from  the  center 
of  the  filament,  we  obtain 

vx(x,0)  a  |  k|  y/ir(x2  +  y2).  (5) 

Notice  that  the  maxiumum  value  is  |ic|/ify,  which  agrees  closely  with  |ic|/2y> 

the  corresponding  velocity  for  a  vortex  ring.  We  may  define  the  mixing  time 

scale  t  ,  as  the  interval  over  which  a  fluid  element  travels  from  (-Sa,  O)  to 
mix 

(Sq,0)  under  the  influence  of  the  residual  flows.  Integrating  Eq.  (5)  gives 

Tmix  =  2*S0<S02/3  +  y2)/|  K|y  •  (6) 

In  Fig.  1,  we  notice  that  the  compression  of  the  hot  region  by  the  flows 

behind  the  shock  has  reduced  the  distance  that  the  fluid  elements  must  travel 
to  reach  the  lower  edge  of  the  heated  region.  Thus  our  tm^x  will  be 
larger  than  the  actual  transit  time. 
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In  the  Appendix,  we  compute  the  integral  in  Eq.  (4)  for  the  case  of  a 
planar  shock  passing  through  a  cylindrically  symmetric  density  depression  like 
that  in  Fig.  1.  We  find  that,  prior  to  the  return  of  the  reflected  shock 
from  the  end  of  the  shock  tube, 

<  “  2  v2(l  -  v2/2W)  S0  toCp^/po)  f 

-  2  v2(W  -  v2/2)  (  t/2)  £n(PM/p0)  f,  (7) 

where  f  <  1  is  a  form  factor  which  varies  with  the  initial  density  profile. 

In  Eq.  (7),  v2  is  the  flow  velocity  behind  rhe  shock  in  the  laboratory  frame; 

W  is  the  shock  velocity;  p^  is  the  ambient  density  of  the  unburned  medium; 

and  Pq  is  the  density  at  the  center  of  the  spherical  flame.  For  a  shock  of 

pressure  ratio  1.3  and  y  -  1.4,  the  Rankine-Hugoniot  relations  give  us  v^3  W/6 

and  W  «  1.12  c  ,  where  c  is  the  ambient  sound  speed  prior  to  passage  of  the 
s  s 

shock.  In  our  simulation  Sn3  2.5  cm,  p  /pna  9.5  ,  and  c  »  3.4  x  104  cm/s. 

u  00  g 

Equation  (7)  gives  us  a  vorticity  strength  of  c<  6  x  lO^cm  2/s.  For  y  <  So, 

the  time  scale  x  .  for  a  fluid  element  to  be  pulled  from  one  side  of  the 
mix 

heated  region  to  the  other  ( Ax  =  2Sq)  under  the  influence  of  the  residual 

vorticity  is  <  900 us.  This  is  consistent  with  the  simulation  results  (Fig.  5) 

and  with  Markstein's  experiment  (Fig.  1),  which  shows  that  fluid  is  pulled 

through  the  center  of  the  burning  zone  prior  to  the  arrival  of  the  reflected 

shock  approximately  600us  after  the  vorticity  is  generated.  We  note  again 

that  in  the  spherical  case  the  fluid  velocity  through  the  center  of  the  vortex 

ring  is  50%  higher  than  that  for  the  two-dimensional  model  which  we  have  used 

here.  Thus  the  shorter  time  scale  (r  .  )  of  Markstein's  experiment  is 

mix 

consistent  with  the  estimate  given  above.  Because  the  region  of  density 
variation  is  no  longer  azimuthally  symmetric,  the  reflected  shock  should 
increase  the  number  of  vortex  centers. 
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NUMERICAL  SIMULATION 


We  have  used  the  two-dimensional  reactive  shock  model  (Oran  et  al. , 

1983),  which  is  based  on  the  code  FAST2D  (Book  et  al . ,  1981),  to  simulate  the 
interaction  between  a  weak  planar  shock  and  a  cylindrically  symmetric  density 
depression  in  air.  The  algorithm  employs  time-step  splitting  in  conjunction 
with  flux-corrected  transport  (FCT)  (Boris  and  Book,  1976)  to  solve  the 
equations  for  conservation  of  mass,  momentum,  and  energy.  The  reactive  shock 
model  also  contains  algorithms  for  chemical  kinetics,  energy  release,  and 
thermal  and  molecular  diffusion.  We  have  calibrated  the  code  extensively 
through  studies  of  shocks  on  wedges  in  nonreactive  fluids  (Book  et  al.,  1981) 
and  studies  of  detonation  cell  structure,  which  required  models  of  energy 
release  (Oran  et  al. ,  1983).  In  the  present  simulations,  we  used  Cartesian 
geometry  and  did  not  include  chemical  reactions  or  diffusive  transport 
processes. 

The  Cartesian  grid  consisted  of  150  x  50  cells  of  dimension  0.2  cm  on 
each  side.  The  grid  remained  fixed  throughout  the  calculation-  Thus  the 
simulated  chamber  was  30  cm  x  10  cm,  as  in  Markstein's  experiment.  We  have 
aligned  the  x-axis  with  the  axis  of  the  shock  tube.  The  simulation  had  inflow 
boundary  conditions  with  v  =  (v^jO)  at  the  side  from  which  the  shock 
propagated,  where  v^  is  the  velocity  behind  the  shock  in  the  laboratory  frame. 
The  other  boundaries  were  reflecting.  A  typical  timestep  was  At  *>  lys, 
where  the  subscript  s  indicates  simulated  time.  Figure  3  is  a  density  contour 
diagram  approximately  90  ys  (simulated  time)  after  the  simulation  began.  The 
density  contours  run  from  3.0  x  10-l+  to  1.26  x  10-3  g/cm3  in  equal  increments. 
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Figure  4  shows  the  initial  density  and  temperature  profiles  more  clearly.  The 
density  profile  had  the  convenient  functional  form 

P(r)  =  p  +  (p  -  p  )/[l+(r/S  )212  .  (8) 

00  0  "  0 

The  shock  pressure  ratio  was  1.3,  again  similar  to  Markstein's  experiment 
(Rudinger,  1958,  and  Markstein,  1964),  and  the  shock  moved  from  top  to  bottom. 
The  ambient  pressure  was  1.0  atm,  and  the  ratio  of  principal  specific  heats 
y(  e  c„/c„)  was  constant  at  the  value  1.4.  The  ambient  density  ( pm)  and 

p  V 

temperature  (TOT)  were  1.17  x  10“ 3  g/cm3  and  300K,  respectively,  while  the 

minimum  density  ( Pg)  and  peak  temperature  (Tg)  were  1.24  x  10~4  g/cm^  and 

2840K.  To  compare  simulation  times  (ts)  with  those  in  Markstein's 

experiment,  one  should  subtract  approximately  150ps  from  the  simulated  time. 

Figures  5-7  show  the  evolution  of  the  density,  pressure,  and  vorticity, 

respectively.  The  x-axis  is  vertical  to  facilitate  comparisons  with  Fig.  1; 

thus  the  shock  initially  propagates  downward  in  the  figures.  There  are 

striking  similarities  to  the  experimental  photographs  (Fig.  1).  Figure  4 

shows  that  the  density  distribution  is  compressed  by  the  shock  and  the  entire 

heated  region  is  pulled  toward  the  bottom  of  the  chamber  by  the  flows  behind 

the  shock.  The  reflected  rarefaction  wave  and  curved  transmitted  shock  are 

clearly  visible  in  Fig.  6.  The  curvature  of  the  latter  occurs  because  of  the 

increased  sound  speed  inside  the  heated  region.  An  interesting  feature  is  the 

Mach  structure  caused  by  reflection  of  the  curved  shock  by  the  chamber.  The 

sequence  of  density  diagrams  clearly  demonstrates  the  effects  of  the  vortex 

pair,  which  has  pulled  ambient  gas  from  the  top  of  the  heated  region  to  the 

bottom  in  a  time  t  .  <  1  ms ,  consistent  with  both  Fig.  1  and  our  theory.  The 

mix  ~ 
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vorticity  contours  in  Fig.  7  show  the  existence  of  two  long-lived  vortices 
which  are  equal  in  magnitude  and  oppositely  directed.  Sometime  after 
the  shock  has  passed  through  the  heated  region  (ts  ~  500ps),  the  simulation 
gives  a  vortex  strength  <  ~  3  x  104  cm2/s,  which  is  consistent  with  the 
theoretical  estimate. 

Another  point  of  considerable  interest  is  the  enhancement  of  the 
vorticity  by  subsequent  shocks.  The  vorticity  has  a  maximum  magnitude  of 
approximately  1.1  x  104  s~  *  just  before  the  reflected  shock  strikes  the  heated 
region  (t  =  722us),  and  1.8xl04  s-1  afterward  (t  =  1.4  ms).  Note  the 
reduction  in  size  and  the  increased  number  of  the  vortex  centers  between  those 
times  in  Fig.  7,  indicating  the  role  of  successive  pressure  waves  in  the 
evolution  of  a  structure  more  like  that  found  in  turbulent  flows. 

Because  the  release  of  chemical  energy  was  not  included  in  the 
calculation  and  because  of  the  resolution  limitations  of  the  calculation, 
several  noteworthy  differences  exist  between  the  experiment  and  the 
simulation.  First  we  notice  in  the  density  contour  diagrams  that  the  heated 
region  spreads  and  that  the  innermost  contour  disappears  by  tg  =  600us, 
indicating  that  cooling  occurs  as  the  ambient  gas  is  pulled  through  the  center 
by  the  vortices.  In  a  reactive  medium,  the  ambient  gas  would  ignite,  at  least 
maintaining  the  temperature  and  producing  pressure  waves,  which  would  hasten 
the  transition  to  turbulence.  The  reactions  would  also  alter  the  reflected 
shock.  Figure  1  (t  =  0.70  ms)  shows  that  the  reflected  shock  wave  actually 
emerges  as  a  dark  band  above  the  vortex  ring.  Finally,  the  limited  resolution 
of  the  calculation  and  the  omission  of  three-dimensional  effects  have 
prevented  much  of  the  cascade  to  small  scale  turbulence,  which  appears  in  Fig. 
1.  Despite  these  limitations,  the  simulation  has  revealed  major  fluid-dynamic 
effects  in  shock  wave  and  flame  interactions. 
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COMPARISON  WITH  RAYLEIGH-TAYLOR  INSTABILITY  THEORY 


Previous  analyses  of  shock  wave  and  flame  interactions  have  used  or  have 
suggested  using  a  modified  version  of  conventional  Rayleigh-Taylor  instability 
theory  (see,  for  example,  Rudinger,  1958,  and  Markstein,  1964).  Indeed, 
Markstein  (private  communication)  has  pointed  out  that  such  an  analysis,  with 
proper  choices  of  parameters,  gives  a  result  similar  to  Eq.  (7).  In  the  usual 
Rayleigh-Taylor  picture,  the  length  scale  of  the  initial  perturbation  is 
infinitesimal  or  "small",  while  the  source  term  in  Eq.  (1), 

S  =  7p  x  7P/o2  (9) 

is  nonzero  for  a  "long"  period  of  time  compared  to  the  growth  rate  of  the 
perturbation.  Markstein  (1964)  modified  this  to  account  for  the  short  time 
scale  of  the  shock  passage  through  a  perturbed  (rippled)  flame  boundary.  He 
demonstrated  the  phenomenon  experimentally  by  passing  a  smooth  flame  through 
a  wire  grid  to  perturb  the  flame  boundary  prior  to  the  arrival  of  a  planar 
shock.  As  expected,  the  small  ripples  in  the  flame  surface  were  amplified  by 
the  perturbed  flow  field  caused  by  the  shock.  These  features  were,  however, 
soon  overwhelmed  by  nonlinear  effects.  Because  the  unperturbed  flame  was 
hemispherical,  the  passage  of  the  planar  shock  produced  large  scale  vorticity, 
which  obliterated  the  smaller  features. 

In  this  paper,  we  have  been  concerned  with  the  more  general,  nonlinear 
description  of  vorticity  generation.  Thus  we  are  addressing  situations  in 
which  the  length  scale  of  the  perturbation  is  finite  (flame  radius)  and  the 
source  term  is  also  nonzero  over  the  time  period  during  which  the  shock  passes 
through  the  entire  flame.  Because  both  scales  are  finite,  we  require  a 
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nonlinear  integration  of  Eq.  (1),  which  covers  the  entire  spectrum  of  inviscid 
vorticity-generating  phenomena.  We  note  that,  even  for  perturbations  of  small 
spatial  and  temporal  extent,  the  analysis  of  the  flow  field  over  long  time 
periods  must  include  nonlinear  effects. 

SUMMARY 


We  have  presented  both  an  analytic  theory  and  numerical  simulations  of 
the  vorticity  generated  by  a  planar  shock  wave  passing  through  a  hot  region 
similar  to  that  of  a  flame.  In  order  to  extract  the  basic  fluid-dynamic 
mechanism  in  this  interaction,  we  have  restricted  our  work  to  two  dimensions 
and  have  excluded  chemical  reactions.  This  is  a  good  approximation  for  the 
calculation  of  the  large  scale  vorticity  distribution,  since  the  laminar  flame 
moves  much  slower  than  the  shock.  Our  analysis  has  shown  that,  after  the 
shock  passes  through  the  flame,  a  residual  rotational  flow  field  remains. 

This  flow  field  is  equivalent  to  that  of  a  vortex  pair  (or  a  vortex  ring  in 
three  dimensions)  characterized  by  a  circulation  <  (Eq.  7).  The  separation 
y  of  the  filaments  is  <  So,  the  radius  of  the  heated  region.  We  have  also 
defined  and  computed  a  mixing  time  scale  for  the  vorticity  distribution  which 
is  consistent  with  Markstein’s  experiment  (Fig.  1)  and  the  simulations. 

Our  numerical  simulations  of  a  planar  shock  interacting  with  a 
cylindrically  symmetric  region  of  reduced  density  used  our  reactive  shock 
model,  which  is  based  on  the  FAST2D  computer  code.  The  computer  model  was 
identical  to  that  used  previously  to  study  detonation  structure  (Oran  et  al., 
1983);  only  the  boundary  and  initial  conditions  were  different.  From  the 
simulations,  we  have  determined  the  time  development  of  the  flow  field  as  well 
as  the  flows  which  remain  after  the  shock  has  passed  through  the  flame.  The 
simulations  agree  well  with  experiment  and  theory,  and  show  that  conventional 
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Rayleigh-Taylor  instability  theory  does  not  provide  an  adequate  description  of 
vorticity  generation.  The  source  of  large  scale  vorticity  for  the  case 
studied  here  is  the  initial  large-scale  misalignment  of  pressure  and  density 
gradients  present  in  the  fluid. 

The  general  framework  presented  in  this  paper  applies  to  vorticity 
generation  by  any  pressure  perturbation  interacting  with  a  density  gradient. 
Examples  of  this  phenomenon  include  sound  waves  or  shocks  passing  through 
flames,  hot  spots,  or  boundary  layers;  the  key  ingredient  is  the  misalignment 
of  the  pressure  and  density  gradients  as  the  interaction  occurs.  We  conclude 
that  pressure  waves  produced  in  one  portion  of  a  reacting  medium  will  generate 
or  enhance  turbulence  as  the  waves  encounter  density  fluctuations  within  the 
medium.  The  scale  lengths  of  the  additional  turbulence  will  be  less  than  the 
spatial  size  of  the  original  density  fluctuations.  We  have  thus  described  and 
demonstrated  a  mechanism  for  producing  turbulence  in  reactive  systems  and  have 
developed  a  framework  for  calculating  the  strength  and  the  scales  on  which  the 
turbulence  is  generated. 
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APPENDIX  -  PLANAR  SHOCK  INTERACTION 


We  will  analyze  the  case  shown  in  Fig.  2  by  considering  (in  two 
dimensions)  a  planar  shock  interacting  with  a  circular  heated  region.  The 
symbol  W  denotes  the  velocity  of  the  shock  through  the  burning  zone,  and  p2, 
v2,  and  P2  are  the  density,  fluid  velocity,  and  pressure  behind  the  shock. 

The  corresponding  quantities  ahf  d  of  the  shock  are  Pp  v^*1  0),  and  Pp  The 

hot  gas  has  a  radius  Sq,  and  the  time  t  =  0  occurs  when  the  shock  reaches  the 

center  of  the  burning  zone.  Here  we  will  assume  that  the  shock  is  weak 

(pressure  ratio  P^P^  *  1*3)  and  that  the  pressure  and  density  ratios  do  not 

change  significantly  when  the  shock  enters  the  burning  zone.  The  circular 
shape  of  the  hot  flame  and  the  elevated  speed  of  sound  in  the  interior  of  the 
flame  cause  the  shock  to  become  curved  upon  passing  through  the  heated  region. 
Our  analysis  will  not  account  for  this  curvature.  We  represent  the  flow  field 
of  the  shock  as  follows: 


p 

°i 

p2  ~ 

pl 

p 

E 

pi 

+ 

P2  ' 

P1 

vi 

_v2  - 

vl_ 

(Al) 


where  s  and  0  define  the  position  of  a  point  relative  to  the  center  of  the 
heated  region,  a  is  a  parameter  defining  shock  thickness,  fa  satisfies  the 
relation 


lim  f 
.  a 
a  -H3 


[  1  -  h(cos  0) 

=  ^h(s  cos9/Wt  -  1) 
tl  -  h(s  cosG/Wt  -  1) 


,  t  =  0 
,  t  <  0 
,  t  >  0, 


(A2) 
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In  this  paper,  we  are  considering  only  variations  in  the  x-y  plane,  so  that 
*  0,  =  5  =  0  and  ?  •  Vv  =  0.  Since  the  integrands  do  not  vary  with  z, 

Eqs.  (1)  and  (A4)  give  us 

dK  2 

~  (t)  /  €  dA(t)  =  /  [( Vp  x  VP)  /p  ]  dA(t).  (A5) 

a  A(t)  Z  A(t)  ~ 

To  determine  the  vorticity  generated  by  a  shock  passing  through  F,  we 
need  to  integrate  eq,  (A5)  over  the  interval  x  during  which  the  shock 
interacts  with  F.  In  the  case  which  we  are  considering,  F  has  a  radius  SQ  ~ 
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2.5  cm.  The  interaction  (shock  transit)  time  will  then  be  quite  short 
compared  to  the  time  scales  of  flame  propagation  and  vortex-induced  rotations  1 
motion.  The  area  A(t)  containing  F  and  nonzero  vorticity  will  change  due  to 
the  motion  of  the  fluid  behind  the  shock  and  the  compression  of  fluid  elements 
by  the  shock.  For  weak  shocks  these  will  also  be  small  effects,  and  we  will 
assume  that  A(t)  =»  Aq,  where  A0  is  the  area  containing  F  prior  to  the  arrival 
of  the  shock.  In  other  words,  we  assume  that  the  Jacobian  matrix  of  the 
transformation  from  dA(t)  to  dA0  is  approximately  the  identity  matrix.  Thus 
we  obtain 


t 

<  (t)  »  /  dt*  /  dx  dy  [(Vpx  VP)  /p2].  (A6) 

Z  .  z 

-  t/2  A 

0 


To  compute  <z(t)  for  a  single  vortex  in  the  upper  half  plane,  we  use 


t  »  » 

<  (t)  =<  /  dt  J  dx  /  dy  [( Vp  x  VP)  /  p2] .  (A7) 

Z  ^  z 


In  Eq.(A7),  we  have  extended  the  limits  of  the  integrals  to  the  entire  half 
plane  and  t  -®  for  convenience  in  performing  the  theoretical  calculations. 
This  is  reasonable  as  long  the  intensity  of  the  interaction  drops  off  rapidly 
as  time  decreases  from  t  =  -t/2. 

We  will  compute  the  vortex  strength  in  Eq.(A7)  in  the  limit  that  the 
shock  has  infinitesimal  thickness  (o  >  0).  We  assume  for  convenience  that  the 
density  profile  of  the  heated  region  is  given  by 
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p(s  ,  t)  =  exp  [-  g(s/S0,  t)  2n(  <>J  Pq)  ] » 


(A8) 


In  computing  <  from  Eqs.  (A7)  and  (4)  (main  text),  we  notice  that  the 
integral  over  the  interval  (-  *>, 0)  for  t  <  0  is  identical  to  that  over  the 
interval  (0,®)  for  t  >  0,  so  that  we  need  only  double  the  result  for  either 
interval  to  determine  the  residual  vortex  strength.  (In  practice  the  time 
integral  is  negligible  outside  the  interval  (-t/2,  t/2),  where  T  is  the  time 
scale  of  the  shock  and  flame  interaction.) 
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Equation  (A7)  now  gives  us 


<  *  2  In  (-A/°dt  /  d  S(  -sin  0)  jdss{-|  6(s-c-g-9  -  1)  x  (All) 


p0  -OO  IT 


2,  ,s  COS  0S 
V2  h(~wi - 1 J  v2s  cose 


With  a  change  of  variables  to  3=3  cos  9/Wt  -  1,  Eq.  (All)  becomes 

2 

pco  r0  *  gg.  v2  h(o)  v2  Wt 

K  ’2  ln  W  idt  i/f(-,1”e>Sls=o  >— (A12) 


Application  of  the  chain  rule  of  differentiation  results  in 

is  =  J*L.  i£  (A13) 

33  „  .  cos  6  3s  u.  » 

3=0  s=Wt/cos9 

which  permits  a  change  of  variables  to  t,  =  Wt/cos0.  Equation  (A12)  is  then 


p  u  v,  h(0)  v,  ” 

(“)  dt  wt  [  is  ~,itd!  H 


=  2  In  ("~)v 2 (W  -  v2/2]  /dt  g(Wt). 
P°  0 


(A14) 


Because  W  >  v2  and  >  p  ,  we  see  that  the  direction  of  the  vorticity  in  Eq. 
(AIM  agrees  with  that  in  Figs.  1  and  5.  We  may  define  a  form  factor  f  by 


f  =  4  /  dt  g(wt), 

T  0 


(A15) 


where  t  =  2  Sn/W,  so  that 


<  a  2  v2  (W  -  v 2/ 2)  (~  )  to  )  f 


(A16) 


-  2  v2  (1  -  v2/2W)  S0  Zn  f. 
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Note  that  we  have  previously  defined  g(s)  =  g(s/S0).  For  a  Gaussian  form, 

g(s/SQ)  =  exp  (-  s2/S02)  (A17) 

<  1 

> 

>  1  (A18) 

f  =  1 .  For  the  Bennett  profile  used  in  the  simulation,  Eq.  (8),  f  will  be 

less  than  /ir/2  (Picone  and  Boris,  1983). 


we  have  f  ®  /ir/2,  while  for  a  square  wave 


g(s/S0)=H 


1  s/S0 
0  s/Sq 
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MARKSTEIN  (1964) 


Fig.  1.  "Interaction  between  a  shock  wave  and  a  flame  of  initially  roughly 
spherical  shape.  Pressure  ratio  of  incident  shock  wave  1.3;  stoichiometric 


n-butane-air  mixture  ignited  at  center  of  combustion  chamber  8.70  ms  before 
origin  of  time  scale."  (Rudinger,  1958  and  Markstein,  1964)  (Reprinted  by 
permission  of  Pergamon  Press,  Inc.). 


THEORETICAL  CALCULATION  NOTATION 


Fig.  2.  Diagram  defining  terms  for  calculation  of  vorutr  strength  in  the 
section  on  theory  and  in  the  appendix. 


NUMERICAL  SIMULATION 

ts-9.4  X  10"®s 
0.0  y(cm)  10.0 


Fig.  3.  Density  contour  diagram  shortly  after  the  simulation  began.  The 
dimensions  of  the  chamber  are  similar  to  those  in  Fig.  1.  Gas  flows  into  the 
chamber  from  the  top  boundary. 
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DENSITY  (10  a  g  cm3) 


NUMERICAL  SIMULATION 


x(cm) 

Fig.  4.  The  density  and  temperature  distributions  at  simulation  time  tg=0 
show  the  shock  at  the  left  and  the  heated  region  in  the  center  of  the  grid. 
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TEMPERATURE  (°K) 


DENSITY  CONTOURS 


0.19  ms  0  38  0  49  0  60 


0.0  (cm)  10.0 


0.72  ms  1.11  1.40 


Fig.  5.  Density  contour  diagrams  at  given  times  tg.  In  the  first  three 
frames  the  six  contours  range  from  3  *  10" 4  to  1.5  x  10“ 3  g/cm3  in  equal 
increments  while  in  the  last  five,  the  contours  range  from  3  x  io-4  to 
1.6  x  10“ 3  g/cm3.  Notice  the  Mach  structure  at  the  channel  boundary  at 


t=384ps. 


In  the  first  three 


Fig.  6.  Pressure  contour  diagrams  at  given  times  t_. 
frames  the  contours  range  from  1.02  x  10 6  to  1.35  x  106  dyne/cm2  in  equal 
increments.  In  the  last  five  the  range  is  1.3  x  io6  to  1.8  x  106  dyne/cm2 
Again  notice  the  Mach  structure  prior  to  reflection  at  the  bottom  boundary 


VORTICITY  CONTOURS 


0.19  ms  0.49 


Fig.  7.  Vorticity  contour  diagrams  at  given  times  tg.  The  ranges  of 
vorticity  values  (respectively  from  earliest  to  latest  times)  are  ±  3  x  103 
s-1,  ±  1.5  x  104  s-1,  -  1.1  x  104s-1,  and  ±  1.8  x  104s_1.  The  number  of 
vortex  centers  and  the  peak  vorticity  increase  after  the  reflected  shock 
passes  through  the  heated  region  (last  frame),  as  anticipated.  Corresponding 
vortices  in  each  frame  are  oppositely  directed  and  oriented  to  pull  fluid  from 
top  to  bottom. 
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